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ABSTRACT 

The sensitivity of the onset and the location of vortex breakdowns in concentrated vor- 
tex cores, and the pronounced tendency of the breakdowns to migrate upstream have been 
characteristic observations of experimental investigations; they have also been features of 
numerical simulations and led to questions about the validity of these simulations. This 
behavior seems to be inconsistent with the strong time-like axial evolution of the flow, as ex- 
pressed explicitly, for example, by the quasi-cylindrical approximate equations for this flow. 
An order-of-magnitude analysis of the equations of motion near breakdown leads to a mod- 
ified set of governing equations, analysis of which demonstrates that the interplay between 
radial inertial, pressure, and viscous forces gives an elliptic character to these concentrated 
swirling flows. Analytical, asymptotic, and numerical solutions of a simplified non-linear 
equation are presented; these qualitatively exhibit the features of vortex onset and location 
noted above. 
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1. Introduction 


The study of the breakdown of the cores of strong longitudinal vortices has been of 
recurring interest for many years. Breakdown may be described as the abrupt change in the 
structure of the core of a vortex, often characterized by the presence of a free stagnation 
point on the axis of the vortex and a corresponding divergence of the stream surfaces near the 
axis, and in some cases a region of reversed flow. In addition to being of theoretical interest, 
vortex breakdown has important technological applications, both aerodynamic and non- 
aerodynamic [1]. In some applications vortex breakdown is beneficial, for example, in leading 
to enhanced mixing or momentum or heat exchange; or destructive, as in degradation of 
aerodynamic performance. Vortex breakdown may play an important role in vortex dynamics 
in general, and perhaps also in transition to turbulence. Given the importance of vortex 
breakdown, it is quite remarkable that the basic underlying mechanism leading to breakdown 
is not yet known or generally accepted! This is not for lack of postulated explanations the 
literature abounds with them. More remarkable yet is that these explanations are so vastly 
different in what they postulate to be the underlying mechanism; these, variously, affirm 
(their principal claimants indicated parenthetically) that: 

(i) vortex breakdown is similar to the separation of a two-dimensional boundary layer 
(Hall [2,3], Mager [4]); 

(ii) vortex breakdown is a consequence of hydrodynamic instability (Ludwieg [5,6]); 

(iii) vortex breakdown is dependent on the existence of a critical state, and is a finite 
transition between states, analogous to a hydraulic jump (Squire [7], Benjamin [8,9]); 

(iv) vortex breakdown is like a solitary wave, or soliton, the result of the trapping of long, 
weakly nonlinear waves propagating in nearly critical swirling flows (Leibovich [10,11]). 

To a greater or lesser degree each of these theories fails to explain fully, or adequately, 
breakdown in all its aspects. Experiments and numerical simulations have failed to resolve 
the matter satisfactorily. How can one account for this failure? A number of reasons present 
themselves. First, experiments and measurements are difficult. There is a tendency for the 
breakdown to migrate back and forth in the test section or on the aerodynamic surface. Mea- 
surements, whether invasive (e.g. hot wires) or non-invasive (e.g. LDV), and flow visualiza- 
tion are difficult to obtain and to interpret in these potentially unsteady, three-dimensional 
flows. Over the years some of the seemingly different postulated physical mechanisms of 
vortex breakdown have been shown to be equivalent or related, leading to the same or sim- 
ilar criteria for breakdown. To the extent that the various theories predict a criterion for 
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vortex breakdown the predictions are quite similar [12], namely that tansup — 1 v/w ~ 57 
(de or v/w xi 1.5, where v and w are the swirl and axial velocities, respectively, so this 
cannot be used to distinguish among or between the various physical mechanisms and their 
corresponding theoretical formulations. 

As for numerical simulations of breakdown using the full Navier-Stokes, or inviscid Euler 
equations, until recently these generally required the assumption of axisymmetry and steadi- 
ness of the flow [13-16]. This is not true of the more recent calculations, for example, those 
by Kuruvila k Salas [17], Spall et al. [18,19], Breuer &; Hanel [20], and Menne k Liu [21], all 
of which solve the unsteady, three-dimensional, fully viscous Navier-Stokes equations (but 
not all of them are necessarily time-accurate solutions). Common characteristics of many of 
the earlier and these most recent numerical simulations, and for that matter of the experi- 
ments as well, are: (i) extreme sensitivity to flow parameters (such as swirl velocity ratio, 
external axial velocity variation, or pressure gradient); (ii) suddenness of breakdown (i.e. no 
evidence of the growth of an instability of the basic swirling flow); and (iii) a tendency for 
the breakdowns to migrate upstream to the initial station or computational cell. In some 
simulations [18] this is avoided by specifying the external axial velocity gradient. 

This last characteristic of these computations is the case, at least, for unbounded or 
“open” flows, e.g. aerodynamic flows and flows in diverging channels. Numerically simulated 
breakdowns, usually in the form of steady axisymmetric bubbles, that occur in bounded or 
confined swirling flows, for example, in finite closed cylinders with one end wall rotating (Lugt 
[22,23], Neitzel [24]) or flow in the gap between rotating spheres (Bar-Yoseph et al. [25]) do 
not exhibit this behavior, nor do experiments on these flows (Escudier [26]). This tendency 
of the breakdowns to migrate to the initial station has caused questions to be raised about 
the validity of such numerical simulations [27]. Another potential concern or question in the 
simulations of breakdown is the appropriate downstream boundary conditions to impose in 
order to “close" the flow domain, a question whether one uses the elliptic Euler or Navier- 
Stokes equations. What makes these issues more troubling is the question of the extent to 
which difficulties with the numerical simulations, including the slow rates of convergence to 
final steady states often encountered, reflect physical reality, or, in other words, are intrinsic 
aspects of the onset and location of breakdown. 

This paper addresses and attempts to answer one of the above questions: the extreme 
sensitivity of the breakdown location to certain flow parameters and the pronounced ten- 
dency for the breakdown to migrate upstream. The theoretical explanation for these effects is 
given in Section 2 where the “quasi-cylindrical equations” (Hall [2]) are analyzed in some de- 
tail. These equations are simplified forms of the Navier-Stokes equations derived by making 
boundary layer-like assumptions and approximations, and like those equations are parabolic, 
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reflecting a strong “time-like” character to the spatial evolution of these thin strongly-swirling 
vortex cores. We show that near breakdown, in the radial momentum equation, which is a 
purely inviscid equation in the quasi-cylindrical approximation, the radial viscous diffusion 
term must be retained. An analysis of this modified set of equations, which would seem, like 
the original equations, to be parabolic since no axial viscous diffusion terms appear, shows 
that in fact they have elliptic-like characteristics. This is explicitly shown by demonstrating 
that the equations may be recast in the form of an elliptic-integro-partial differential equa- 
tion for the axial velocity variation. In Section 3 a simplified, but still non-linear, equation 
is employed as a model of the full integro-differential equation and its full solution given. 
Included in this section is also a matched asymptotic expansion of the equation in the limit 
of large Reynolds numbers. (This solution yields simpler expressions for the axial velocity 
variation near breakdown and its explicit dependence on flow and fluid parameters.) This 
section concludes with an order-of-magnitude analysis which explicitly relates the constant 
appearing in the simplified equation to the primary flow parameters: the core Reynolds 
number, the flow divergence, the swirl velocity ratio, and the (non-dimensional) core circu- 
lation. In Section 4 explicit results of the exact analytic solution presented in Section 3 for 
a representative range of values of the flow parameters are presented and discussed. 


2. Mathematical Model 


We assume the flow is steady, axisymmetric and laminar. The governing Navier Stokes 
equations in a cylindrical coordinate system r, 0, z, with corresponding velocity components 
u, V, w (the overbars denote dimensional variables) are 


1 d . . dw 

Fa? (r “ )+ a? = 0 ' 


_du _du v 2 1 dp \ d (l d 

U dr W dz r pdr + V dr \r dr 


(r u) + 


d 2 u' 


_dv dv u v 

u — + w— + — = v 
or oz r 


d (ld_ \ dhv' 

dr ( r dr r *7 + dz 2 


_dw _dw 1 dp 

u dF + w d¥ = ~di + '' 


'1 d (_ dw\ , 
rdr V dr ) 


d 2 w' 
dz 2 


( 2 . 1 ) 

( 2 . 2 ) 

(2.3) 

(2.4) 
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We define the following non-dimensional variables 


f L\ _U_ _ p - Poo 

U/wy p ~ P ws ' 



(2.5) 


where L is a characteristic scale in the z-direction, 6 is a characteristic core radius, W 0 a 
characteristic axial flow velocity (e.g., free-stream velocity at the edge of the core), the 
uniform static pressure far from the vortex, p the density, and v the kinematic viscosity. In 
terms of these quantities the governing equations become 


\_d_ 

r dr 


/ x dw 
(ru) + = 0, 


(2.6) 


/ 8 \ 2 [ du du 

UJ [ U Tr +W Yz 


V 

r 


dp 1 
dr + Res 



d_ 

dr 






(2.7) 

(2.8) 

(2.9) 


where Res = WqS/v is the Reynolds number based on the characteristic core radius. 

Equations (2.6) - (2.9) contain two parameters: 8jL and Re$. At large distances upstream 
of breakdown the vortex develops on a length scale L which is much larger than 8 , so 8/ L <C 1 
in this region. The core Reynolds numbers for which vortex breakdowns occur are large, 
but only modestly so, of the order of hundreds ([14-16], [18]). We see immediately that 
if 8/L = t <C 1 that the axial viscous diffusion terms in the radial, azimuthal, and axial 
momentum equations, (2.7) - (2.9), are negligible compared to the other terms in these 
equations. If we take the formal limit e — + 0, with Re$ held fixed, and without assuming 
that eRes is either large or small, then (2.6) - (2.9) reduce to 


1 d . . dw 

rYr [Tu) + ZZ = °’ 


( 2 . 10 ) 


V 


2 


r 


dp 


dv dv uv 

Uq- + W— + = 

or oz r 



( 2 . 11 ) 

( 2 . 12 ) 
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(2.13) 


dw dw 

u ^~ + w ~5~ 
dr dz 


dp . 1 

dz Re s 



Equations (2.10) - (2.13) are the so-called quasi-cylindrical equations, derived and extensively 
analyzed by Hall [3]. In contrast to Eqs. (2.6) - (2.9) which are elliptic, these quasi-cylindrical 
equations are parabolic equations, which can be solved by time-marching schemes. They 
differ from the ordinary boundary layer equations only by the inviscid radial momentum 
equation, (2.11), which contributes a non-zero normal pressure gradient. These equations 
were solved numerically by Hall [27] to predict the location of the onset of breakdown. 
Because upstream influence has been neglected they cannot be expected to predict the details 
of the breakdown itself. 

As we approach breakdown the scale of z-derivatives becomes comparable to that of r- 
derivatives, so 8/L, even if small, is not expected to be infinitesimal. Hence, returning to 
Eqs. (2.6) - (2.9), for this near-breakdown region we again drop all terms of 0[(8/ T) 2 ], but 
keep in addition to the (L/8)ReJ l terms also terms of 0(8/ L). (In other words, we are 
assuming that 6/L is small, but not infinitesimal, so terms of 0[(8/ L) 2 } may be neglected, 
but not those of 0(8/ L).) There is only one such term in Eqs. (2.6) - (2.9), the second term 
on the right hand side of (2.7); this equation then becomes 


2 dp 1 (6\ d ( Id, A 


(2.14) 


The other equations remain unchanged, viz. (2.10), (2.12) and (2.13). (Eq. (2.14) as a 
replacement for (2.11) in the near-breakdown region may be more formally justified by an 
appropriate rescaling of variables or by a multiple-scale argument.) In dimensional form the 
new set of reduced equations for the near-breakdown region is then 


Id dw 

^ (rtt)+ a7 = 0 ' 

(2.15) 

du du v 2 l dp d ( 1 d \ 

u Tr + w Tz - 7 = —pFr + V F (y* (ru) j ' 

(2.16) 

dv dv uv d (\ 1 \ 

U d? + W rz + V = ,, Vr{rFr {rv) )' 

(2.17) 

dw dw l dp v d ( dw\ 

U !fr + W lh ~ ~~pdz + V -fc ) ' 

(2.18) 


These equations differ from the full axisymmetric Navier-Stokes equations only by the omis- 
sion of all the axial viscous diffusion terms and the inertia terms in the radial momentum 
equation. They differ from the quasi-cylindrical equations, (2.10) - (2.13), only by the inclu- 
sion of the radial viscous diffusion term in the radial momentum equation. Since, as for this 
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modified set of reduced equations, no axial viscous diffusion terms appear in (2.15) - (2.18) 
we would expect this set to be parabolic as well. As we shall see below, however, this is not 


so. 


Combining (2.15) and (2.16) we can write 

,2 


v 1 dp d 2 w 

t p dr drdz 


(2.19) 


If we integrate this with respect to r from 0 to r, and then differentiate the resulting expres- 
sion with respect to z, it becomes 

d 2 w , d 2 w . 


dp 

d~z 


) „“(&) “"'i vk dr+ 




r=0 


dz 2 ' dz 2 

Multiplication of (2.17) by 2u, and combining terms, leads to the expression 


9 , 2 2 x , d(r 2 v 2 ) 2 

u—(r v ) + w ~ 


d 


= 2 vr v — 
or 


\ 9 ‘ 
r¥r {rv) 


( 2 . 20 ) 


( 2 . 21 ) 


dr v ' ' ' ~ dz 

which can be used to eliminate dv 2 /dz in the integral term in (2.20), resulting in, dropping 


/ v \ 9 

Id, . 


'd 2 w 

92 W j 

(w) d~r 

~jr( rv ) 
r or 

dr+fi 

dz 2 

dz 2 |r=0 . 


, f 1 ( “) #(rVMr-2„ ['-(-) l 
dz dz ^ Jo r 3 \wj dr ^ Jo r \wJ dr 

( 2 . 22 ) 

Equation (2.22) exhibits the various mechanisms that can lead to the development of an 
adverse pressure gradient along the core axis and consequent axial retardation of the flow. 
The first term on the right hand side is an inviscid term and was first identified in this context 
by Hall [3] from an inviscid analysis of the quasi-cylindrical equations. It is in the nature of a 
flow divergence-swirl interaction term, demonstrating how, for u > 0, arising from whatever 
cause, and positive dk/dr, where k — rv , which is generally the case, the pressure gradient 
along the axis can exceed that for larger r, and in particular that at the core edge. If T and 
R are characteristic measures of the core circulation (i.e., rv evaluated at the core edge) and 
radius, and we set u/w = cv, a measure of flow divergence, an order-of-magnitude analysis 
of this integral, evaluated at the core edge, yields (Hall [28]) 

dp dp. paT 2 


dz W =° dz' T=n R 3 ' 

Thus, if a > 0, i.e. the stream surfaces diverge in the core, then the pressure gradient on the 
axis exceeds that along the core edge by an amount that can be very large for a concentrated, 
intense vortex core. 

The two new terms on the right hand side of (2.22) are both viscous terms. The first of 
these, the second one on the right hand side, represents the interaction of swirl and radial 


dp 


(2.23) 
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diffusion of swirl momentum, and also acts to raise the pressure on the axis. The last 
term on the right has the form of an axial diffusion term, which has arisen not from the 
inclusion of axial diffusion terms in the governing equations, because there are conspicuously 
absent from Eqs. (2.15) - (2.18), but from the continuity equation. Thus while the original 
quasi-cylindrical equations are parabolic, their modification, by the inclusion of the radial 
diffusion term in the radial momentum equation, has the effect of introducing a term which 
has the appearance of an axial diffusion term! This term gives an elliptic character to 
the equations. More precisely, (2.22) is an elliptic-integro-partial differential equation. In 
addition to showing how upstream influence enters the vortex breakdown problem, Eq. (2.22) 
also exhibits explicitly how flow divergence, a = u/w , and the swirl ratio, v/w, independently 
affect the axial pressure gradient. 

We can see more immediately the effect of the terms on the right hand side of (2.22) on 
the axial velocity by using (2.18) to eliminate dpjdz from (2.22), obtaining 


dw 

‘"aT 1 "-’ 


’ld_ 

r dr 



T = 0 


dw dw\ 

u a7 + “’aT j 


v d f dw' 
r dr 1 dr 


= - Jo h (z) Tr {rW)d: 


(2.24) 


. r \ ( v\ d \\ d . . 
+2v - - — --r-(rv) 

Jo r \wj dr r dr 


dr — v 


d 2 w d 2 


w, 


dz 2 dz 2 


r=0 


Equation (2.24) contains as unknowns the three velocity components. To complete the 
formulation of the problem we introduce the stream function 0, defined such that 


1 dip 
U = —~r~dz' 


1 dxp 
w = — 
r dr 


(2.25) 


which makes the continuity equation automatically satisfied. Substitution of (2.25) into the 
azimuthal momentum equation, (2.17), results in 


, dv , , dv d(\d \ 


(2.26) 


Eliminating u and w in (2.24) using (2.25) yields 
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Equations (2.26) and (2.27) are the governing equations for the unknowns j/> and v subject 
to appropriate boundary conditions. We note again that (2.27) is an elliptic-integro-partial 
differential equation. 


3. Simplified Model 


We noted earlier, in a discussion of expression (2.23), that if the flow divergence a = 
u/w > 0, then it follows that (dp/dz) r= o > ( dp/dz) T= R , i.e., flow divergence in such vortex 
flows can lead to a less favorable pressure gradient along the axis compared to that at the 
edge of the core. Equation (2.23) is based on a purely inviscid argument, nor is anything 
said about the origin of the streamline divergence. Within the inviscid assumption this may 
most readily be assumed to arise from an imposed adverse pressure gradient at the core 
edge, leading to axial flow retardation, and from continuity to u > 0. Could viscous effects 
either cause or alter the degree of axial retardation, independent of purely inviscid effects. 
Grabowski and Berger [14] conjectured that this could come about as follows. Viscous 
dissipation of swirl would lead to a decrease in v (Eq. (2.17)), so dv/dz < 0. It then follows 
from (2.11) that there will be greater adverse axial pressure gradient and axial velocity 
retardation near the axis than at the core edge. This, according to continuity, (2.15), will 
lead to a significant amount of outflow (u > 0). This outflow, in turn, because of conservation 
of angular momentum, will lead to a decrease in v (because r increases). The combination of 
outflow, and therefore flow divergence, and diffusion of swirl, leading to a decrease in v and 
hence making dvjdz < 0 will combine to make for a large adverse axial pressure gradient 
and resulting axial retardation near the axis. This picture of the flow closely parallels that 
of Hall [28]. 

The analysis presented in Section 2 attempts to capture this physics. The resulting for- 
mulation is, however, too complicated for the essentially new elements, those that distinguish 
this formulation from the quasi-cylindrical model, to be evident. In an attempt to illustrate 
these, we consider in this section a simplified version of the analysis of Section 2. 

If we evaluate (2.24) at r = R, the edge of the core, noting that the viscous terms on the 
left hand side will be small there, as well as u dw/dr , w dwjdz , and d 2 wjdz 2 , it reduces to 


w 


dw . 

lh' 


- = -/“i(")r< A > +2 '/ S1 (-) 

Jo \wj or Jo r \w J 


d_ 

dr 


Id 

-7r(™) 

r or 


d 2 w. , . 

dr + v——\ T = 0 . (3.1) 
dz z 


We anticipate that in situations where breakdown occurs the sum of the two integrals on the 
right hand side of (3.1) is negative since this will lead to retardation on the axis (dw/dz < 0). 
Because these are integral terms they are less sensitive to the details of the velocity profiles, 
and as our greatest simplification we set the sum of these terms to a (negative) constant 
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— A, say. Equation (3.1) then becomes 


or, for simplicity 


d 2 w. dw. , . n , 

V-—\ T= o - W~~ r =0 & A, {A > °) 

<9z 2 


d 2 ic dw 

v-nr - w— = > 1 , 


dz 2 dz 

where w in (3.3) is understood to be tc| r= o, the axial velocity on the axis. 
First, we note that for small or negligible viscosity, 

dw 


(3.2) 


(3.3) 


with solution 


1 


W— Ri —A, 

dz 


w 2 ~ — Az + const. 


(3.4) 


If w = W c at z = 0, the uniform axial core velocity far upstream, then the const, in the 
above is / 2 and 


w » W c 1 - 


2A 

w~r 


1/2 


For small z this gives 


w ~ W c[ l -wl z + 


(3.5) 


(3.6) 


a linear decrease of axial velocity on the axis. This is in qualitative and quantitative agree- 
ment with Hall’s calculations [28] using the quasi-cylindrical equations. 

3.1. Exact solution of non-linear model equation 

Here we present the exact solution of non-linear Eq. (3.3). Integration yields immediately 

— uf 2 = f — 'j 2 + a (a = const.). (3.7) 

dz 2v \v J 

This is a Riccati equation. It can be transformed to a second-order linear ordinary differential 
equation by setting 


Equation (3.7) becomes 


d 2 y 1 (A 


y'(z) 

~ 2u ( i - 

y \ z ) 

(3.8) 

A \ 

—z + a)y = 0, 
v ) 

(3.9) 


which is an Airy equation. It can be put in standard Airy equation form by introducing the 
new independent variable 

y/ 3 ,M ' 


v 


whereupon (3.9) becomes 


(3.11) 


The general solution of (3.11) is 


d 2y . n 

dP~ Zy = 0 ' 


y(z) = CiAi(z) + c 2 Bi(z), 


(3.12) 


where C\ and C 2 are arbitrary constants and Ai and Bi the Airy functions. In terms of the 
original dependent variable w the solution of (3.3) is then 


w(z) = (4 Av) 1 ! 3 


kAi'(z) + Bi’(z) 
kAi(z) + Bi(z ) ’ 


(3.13) 


where the arbitrary constants for the general solution of the second-order differential equation 
are a and k. 

As boundary or initial conditions for (3.13) we have 


at z = 0. (3.14a, b ) 

The first of these has already been discussed. As for the second, far upstream, at z — 0, 
the change in w with z may be small in the absence of an imposed pressure gradient, but 
this would not be so if this was not the case. Thus, for example, by choosing ( dw/dz)o < 0 
we can consider cases where an adverse pressure gradient, say as resulting from a diverging 
channel, is imposed at the core edge. 

Applying first (3.14b) to (3.7), we find immediately 


w = W c 


dw ( dw\ 

Tz^ (dj) ,glven 



(3.15) 

Since, when z = 0, 

( v y /3 

*= \ 2A 2 ) « = - 0 , say, 

(3.16) 

so, applying boundary condition (3.14a) 


kAi(z 0 ) + Bi(z 0 ) 

(3.17) 

Solving for k, 

W c Bi(z 0 ) - (4Ai>y/ 3 Bi'(z 0 ) 
W c Ai{z Q ) - {AAiiyi 3 Ai'{z 0 ) ' 

(3.18) 
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This completes the determination of the two constants a and k appearing in the general 
solution, (3.13), of the non-linear equation (3.3). The solution (3.13) depends on the four 
(dimensional) parameters A,i/,W C , and (dw/dz)o- 
We note that if (dw/dz ) o = 0, then 

a = -^Wl (3.19) 

and the solution simplifies accordingly. 


3.2. Asymptotic solution for Re$ — » oo 


We consider here, using matched asymptotic expansions, the asymptotic solution of the 
non-linear ordinary differential equation (3.3) for vanishingly small viscosity u, i.e., in the 
limit Res — * oo. 

We begin with the outer expansion, which follows immediately from the first integral of 
(3.3), namely the expression (3.7), dropping, to lowest order, the first term on the left, so 


— w 2 ~ — Az + const. 

2 


(3.20) 


If w = W c at z = 0, then the const. = (1/2)W C 2 , and the leading order outer solution is 

1/2 

| , (3.21) 


Tir / 2 A 

W ~ W r I 1 — TT7TZ 


W? 


which is exactly the earlier (3.5). 

To determine the scale of the inner variables, z and w, we set 


w = Dv a w (a < 0) 


where we anticipate both a,/? < 0. The differential equation (3.3) becomes 

u W-«+l^L __ D -l u -2a+P w ^_ = AD 

dz* <iz 

The two terms on the left hand side will be of the same order in v if 


(3.22a, b ) 


(3.23) 


ft + a = - 1 . 

If we choose a and (3 such that the single term on the right hand side vanishes in the limit 
as v — i 0, then the resulting inner solution cannot be matched with the outer solution. 
Therefore we must choose 

—2 a + (3 — 0, or (3 = 2a. 
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It then follows that 

a = -1/3, fi = - 2/3. 
Thus the scaling of the inner variables is 

Dw 


(3.24) 


w 


z = 


1/1/2’ 

(W? 


2A 


or w = D V/ 3 tn 


— z') f/ 2/ ^ 3 , or z = — j 

/ 2A 


(3.25) 


Equation (3.23) becomes 


cPw _._i _dw 

— + d ^- = ad. 


(3.26) 


This equation for the inner region is identical, apart from a change of variables, with the 
full equation, (3.3). The following analysis therefore proceeds along the lines of that given 
in Section 3.1. 

The first integral of (3.26) is 


die 1 , , N 

— - -| w = ADz + a, (a = const.) 

dz 2D v 


(3.27) 


This Ricatti equation can be converted to a second order linear equation by the substitution 
W oc y' jy. The solution of the resulting Airy equation is 


where 


n kAi’(z) + Bf(S) 

w{z) ~ (4i4) D kW) + W) ’ 


(2A 2 )" 1 * 3 

B = i£±l —(ADz + a), 


(3.28a) 


(3.286) 


a and k being arbitrary constants, to be determined by matching of the inner and outer 
solutions. 

To match the outer solution (3.21) and the inner solution (3.28), we begin by rewriting 
the outer solution (3.21) in inner variables and expanding for v — ► 0, with z held fixed, 
obtaining 


w, 


outer 


2AzS'\ 


Rewritten in outer variables the inner variable z is 


z = 


^-2/3 

2 J / 3 D 


AD {^- z )^' 3+a 


(3.29) 


(3.30) 
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For reasons discussed earlier A > 0, and we can assume that D , a scale factor for w, is also 
positive. It then follows that in the limit // — » 0, with 2 held fixed, that 


W 2 

z — > 00 if — 7 — z > 0 , 
2A 


W 2 

z — > —00 if — 7 — z < 0 . 

2A 


(3.31) 


If the second of these alternatives is the case, Ai(z) and Bi(z ) oscillate and the inner solution 
(3.28) cannot be matched to the outer solution, (3.29). Therefore we assume the first case, 
that z — * 00 as v — > 0. For asymptotically large positive arguments Ai and Bi have the 
following asymptotic expansions 


as z — > 00 < 


m(z) ~ 

A>V) ~ 

r 


[ Bi(z) ~ 4= i' 1 / 4 e 2i3/2/3 , 

1 

Bi'(z) ~ -^z 1 / 4 e 2 | 3 / 2 /3. 
v 7r 


(3.32) 


It then follows that in this (outer) limit that (3.28a) behaves asymptotically as 


w{z) ~ (4 AY /3 Dz l/2 , 


(3.33) 


so, according to (3.25) and (3.28b), 

( 4/1) 1 /V/ 3 z 1 / 2 

(2y4 ) -1 / 3 


w i, 


( 4 / 1 ) 1 / 3 z / 1 / 3 


1 1/2 


or. 


D 

2a\ 


(ADz + a) 


^ inner ~ ( 2 ^ + jj) ’ /2/3 ‘ ( 3 - 34 ) 

The variable part of the inner and outer solutions, (3.29) and (3.34) agree exactly. 

3.3. Exact solution in the limit of vanishingly small viscosity 

We consider the limit of the exact solution given in Section 3.1 as v —* 0. In this limit, 
from (3.15) and (3.16), 

a -* -^-W 2 -> 00 , (3.35a) 


i w ? - 00. 

V \ 1 / 3 / 1 




(3.356) 
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From (3.18), using the asymptotic expressions (3.32), 


lim k = 0. 

i/-» o 


(3.35c) 


From (3.10) 


Ura.-oi - -(© ' (r-^) 
-» -(2TV)- 1/3 (W 2 - 22b) 


(3.36a) 


' +oo, if W c 2 > 2Az 

-» i (3.366) 

k -oo, if W 2 < 2Az. 

We assume that IF 2 > 2Az, so z —> +oo as v — > 0. It then follows from (3.13), using (3.35c), 
(3.32) and (3.36a), that 


lim w(z) (4/b) 1/3 i 1/2 = (42b) 1/3 ^V)' 1 / 3 (w 2 - 2Az) 

-+ {W?-2Azy! 2 , 


1 1/2 


(3.37) 


which is identical with the outer solution, (3.21), presented in Section 3.2 as part of the 
asymptotic solution for Res — » oo, and with the limiting small viscosity solution given 
earlier (Eq. 3.5). 


3.4. Non-dimensionalization 


We begin our non-dimensionalization of (3.3) by first considering the quantity on the 
right hand side, A, which stands for the (negative of the) two integral terms on the right 
hand side of (3.1). We consider each of these in turn. The order of magnitude of the first of 
these can be estimated as follows: 




(3.38) 


where rv = k ~ T, the characteristic circulation of the vortex, 8 the characteristic core 
radius (~ R), and u/w ~ a, a characteristic flow divergence. 

The magnitude of the second integral on the right hand side of (3.1) is similarly estimated: 


2v 



’\d ' 


dr ~ 2 " (£) © 8 = 2 " (JO r> 


(3.39) 


where v/w ~ s, a characteristic swirl velocity ratio. (Both integrals have dimensions of 
LT~ 2 .) Unlike (3.38) which is generally negative as the vortex approaches breakdown, the 
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sign of the second integral is not obvious. To allow us to consider the effect of the integral 
being negative (driving the vortex to breakdown) or positive, we shall attach to it an arbitrary 
sign, denoted by sgn. Then 

A~^r 2 + 2 sgn 1/ (£) r, (3.40) 


and the differential equation (3.3) becomes 


cPw 


V 


dz 2 



( i ) r ' 2 + 2 sgn 



(3.41) 


We now non-dimensionalize (3.41) by defining non-dimensional variables and parameters 


_ z w 

J= l' W ~ We’ 


7= jl Res = m 


(3.42) 


SW r 


with which definitions (3.41) becomes 


1 cPw dw 2 2 

Re s dz 2 W dz ~ ^ + Re s 


sgnsy. 


(.3.43) 


The solution of (3.43) depends on the four parameters Reg, the core Reynolds number, 
a(~ u/u;), the flow divergence, s(~ vj w), the swirl velocity ratio, and 7(= T/^Wc), the 
non-dimensional flow circulation. 


4. Calculations and Results 

We have calculated a number of solutions of Eq. (3.43), by evaluation of the exact 
solution given in Section 3.1, in order to exhibit the dependence of the solutions on the 
four parameters Reg, a, s , and 7 which appear in the equation, and the value of the initial 
velocity gradient ( dw/dz) 0 ■ Figures 1 and 2 show the non-dimensional axial velocity variation 
w(z) with Reg (for convenience the bars over the non-dimensional variables introduced in 
Section 3.4 have been removed) for the same values of (dw/dz)o(= —0.1) and c*(= 0.005) 
but different values of s(= 0.5 and 1.0) and 7(= 1.0 and 1.5). (We anticipate that negative 
values of (dw/dz)o, corresponding to an adverse axial pressure gradient, are most likely to 
lead to breakdown; this is borne out by sample calculations of such cases. Also, in all the 
calculations reported here the sgn in (3.40) is taken to be positive, implying that the integral 
in (3.39) is negative, as is that in (3.38), increasing the tendency of the vortex to undergo 
breakdown.) We note that for both cases in Figs. 1 and 2 the vortex undergoes breakdown, 
i.e., the axial velocity becomes zero, for all Res in the range 50 < Re$ < 1000, except 
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possibly for Re s = 1000. Moreover, the point on the axis where stagnation or breakdown 
occurs moves rapidly upstream as Res increases, reaching asymptotically a position close to 
the initial station. Finally we note that the curves in Figs. 1 and 2 are remarkably similar for 
the same Re s , suggesting that the results are insensitive to the values of s and 7. Calculations 
confirm this for other values of the parameters as well. Figure 3 shows the axial velocity 
for the same Res range and the same s and 7 as in Fig. 2 but a smaller value of a and a 
larger negative initial axial velocity gradient. We note in comparing Fig. 3 to 2 the much 
more rapid decrease in axial velocity as Re , 5 increases. This is primarily a consequence of 
the much larger negative value of (dw/dz ) 0 in Fig. 3. (That the change in a has only a small 
effect can be seen, for example, from Fig. 4 for the same values of all the parameters as in 
Fig. 3 except for ct which is ten times as large, 0.01 as compared to 0.001. There is very 
little discernible difference between the curves in Figs. 3 and 4.) Figure 5 shows the same 
curves for the same values of the parameters as in Fig. 2 except that for Fig. 5 \(dw / dz)o\ is 
a tenth as large. The consequence of this is that compared to Fig. 2 the curves are displaced 
further downstream, and even more significantly, for the smaller Res the axial flow does not 
stagnate but remains relatively constant and then increases with axial distance! Again this 
illustrates the importance of the value of the initial axial velocity gradient. 

Because of the particular significance of the value of the axial location where the axial 
velocity becomes zero, denoting the location of breakdown, we have used the method of con- 
tinuation [29] to determine the variation of this point as a function of the various parameters. 
The results above, as well as others unreported here, demonstrate that the point at which 
the axial velocity becomes zero, denoted by z v . h ., varies primarily with the parameters Re s 
and (dw/dz)o- Figure 6 shows the variation of z v _b. with Res for four values of (^dxv j dz) 0, 

the other parameters have the values a = 0.01, 7=1, and s = 0.5. We see that for all these 
values of the initial gradient the location of axial stagnation moves rapidly upstream, and 
that while the exact values of z v _b. asymptote for larger Res there are significant differences 
for smaller Res - Figure 7 shows, for Re = 500, a = 0.015, 7=1, and s = 0.5, the variation 
of z v ,b. with —(dw/dz) 0 . We note the very rapid upstream migration of the location of axial 
breakdown as the initial axial velocity gradient becomes more negative, i.e., larger initial 
adverse pressure gradient. 

5. Conclusions 

The principal aims of this paper were to show the origin of the pronounced tendency of the 
location of vortex breakdown to migrate upstream, an effect characteristic of most numerical 
simulations and noted in some experiments. Based on an order-of-magnitude analysis of the 
viscous Navier-Stokes equations a modified set of governing differential equations for swirling 
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flows has been derived. A theoretical analysis of these equations shows that unlike the 
parabolic quasi-cylindrical equations, which they resemble, these equations have an elliptic 
character. This analysis also exhibits how the location of axial stagnation, interpreted here 
as a signal of breakdown, depends explicitly on a number of parameters: the core Reynolds 
number, the flow divergence, the swirl velocity ratio, and the strength of the vortex. A 
simplified model for the axial velocity variation is presented and solved exactly. Numerical 
results obtained from this solution for a range of the governing parameters show that there 
is a very pronounced tendency of the breakdown location to move upstream as the core 
Reynolds number increases, or the initial adverse pressure gradient increases, in qualitative 
accord with most numerical and experimental simulations. 
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Figure 1. Variation of axial velocity with axial distance for different core Reynolds numbers; 
parameters: a = 0.005, 7 = 1.5, s — 1.0; initial axial velocity gradient, ( dw/dz)o = —0.1. 



Z 

Figure 2. Variation of axial velocity with axial distance for different core Reynolds numbers; 
parameters: a — 0.005, 7 = 1, s = 0.5; initial axial velocity gradient, (div/dz ) 0 = —0.1. 







z 

Figure 3. Variation of axial velocity with axial distance for different core Reynolds numbers; 
parameters: a = 0.001, 7=1, s = 0.5; initial axial velocity gradient, (dw/dz) 0 = -0.15. 



Z 

Figure 4. Variation of axial velocity with axial distance for different core Reynolds numbers; 
parameters: a = 0.01, 7=1, s = 0.5; initial axial velocity gradient, (dw/dz) o = —0.15. 
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Figure 5. Variation of axial velocity with axial distance for different core Reynolds numbers; 
parameters: a = 0.005, 7 = 1, s = 0.5; initial axial velocity gradient, ( dw/dz) 0 = —0.01. 
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Re 6 

Figure 6. Location of the axial vortex breakdown point as a function of core Reynolds 
number for different initial axial velocity gradients; parameters: a = 0.01, 7=1, 3 = 0.5. 
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